library(fitdistrplus)
Cargando paquete requerido: MASS
Cargando paquete requerido: survival
library(tidyverse)
── Attaching core tidyverse packages ────────────── tidyverse 2.0.0 ──
✔ dplyr 1.1.4 ✔ readr 2.1.5
✔ forcats 1.0.0 ✔ stringr 1.5.1
✔ ggplot2 3.5.1 ✔ tibble 3.2.1
✔ lubridate 1.9.4 ✔ tidyr 1.3.1
✔ purrr 1.0.4 ── Conflicts ──────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag() masks stats::lag()
✖ dplyr::select() masks MASS::select()
ℹ Use the ]8;;http://conflicted.r-lib.org/conflicted package]8;; to force all conflicts to become errors
library(queueing)
library(kableExtra)
Adjuntando el paquete: ‘kableExtra’
The following object is masked from ‘package:dplyr’:
group_rows
dfRaw <- vroom::vroom(file = '../RawData/data_sst_reto_coppel.csv', delim = ',')
New names:Rows: 18510772 Columns: 10── Column specification ──────────────────────────────────────────────
Delimiter: ","
chr (5): Segmento, caja, tienda, status, estado
dbl (4): ...1, hora_llegada, hora_llamado, hora_salida
date (1): Fecha
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
str(dfRaw)
spc_tbl_ [18,510,772 × 10] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
$ ...1 : num [1:18510772] 821243 821244 821245 821246 821247 ...
$ Fecha : Date[1:18510772], format: "2023-01-01" ...
$ Segmento : chr [1:18510772] "retail" "retail" "retail" "retail" ...
$ hora_llegada: num [1:18510772] 0.503 0.504 0.504 0.504 0.504 ...
$ hora_llamado: num [1:18510772] 0.503 0.504 0.504 0.505 0.508 ...
$ hora_salida : num [1:18510772] 0.504 0.511 0.505 0.509 0.51 ...
$ caja : chr [1:18510772] "caja_a" "caja_a" "caja_b" "caja_b" ...
$ tienda : chr [1:18510772] "Tienda_A" "Tienda_A" "Tienda_A" "Tienda_A" ...
$ status : chr [1:18510772] "Atendido" "Atendido" "Atendido" "Atendido" ...
$ estado : chr [1:18510772] "Ciudad de México" "Ciudad de México" "Ciudad de México" "Ciudad de México" ...
- attr(*, "spec")=
.. cols(
.. ...1 = col_double(),
.. Fecha = col_date(format = ""),
.. Segmento = col_character(),
.. hora_llegada = col_double(),
.. hora_llamado = col_double(),
.. hora_salida = col_double(),
.. caja = col_character(),
.. tienda = col_character(),
.. status = col_character(),
.. estado = col_character()
.. )
- attr(*, "problems")=<externalptr>
# cargar funciones auxiliares
source('../Code/R/utils.R')
Veamos como estan los datos
glimpse(dfRaw)
Rows: 18,510,772
Columns: 10
$ ...1 <dbl> 821243, 821244, 821245, 821246, 821247, 821248, 821249, 82125…
$ Fecha <date> 2023-01-01, 2023-01-01, 2023-01-01, 2023-01-01, 2023-01-01, …
$ Segmento <chr> "retail", "retail", "retail", "retail", "retail", "retail", "…
$ hora_llegada <dbl> 0.5034653, 0.5035810, 0.5037894, 0.5038472, 0.5038935, 0.5052…
$ hora_llamado <dbl> 0.5034769, 0.5036157, 0.5038009, 0.5051319, 0.5080718, 0.5093…
$ hora_salida <dbl> 0.5035810, 0.5111736, 0.5051204, 0.5093333, 0.5099815, 0.5099…
$ caja <chr> "caja_a", "caja_a", "caja_b", "caja_b", "caja_c", "caja_b", "…
$ tienda <chr> "Tienda_A", "Tienda_A", "Tienda_A", "Tienda_A", "Tienda_A", "…
$ status <chr> "Atendido", "Atendido", "Atendido", "Atendido", "Atendido", "…
$ estado <chr> "Ciudad de México", "Ciudad de México", "Ciudad de México", "…
los tipos de datos son correctos, aunque las horas de llegada, llamado y salida vienen en double
Falta transformar a la hora en sitema sexagesimal para una mejor apreciacion de los datos
dfRaw |>
summarise(across(everything(), \(x)sum(is.na(x)))) |>
kbl() |>
kable_classic()
| ...1 | Fecha | Segmento | hora_llegada | hora_llamado | hora_salida | caja | tienda | status | estado |
|---|---|---|---|---|---|---|---|---|---|
| 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
veo que no hay ningun NA, y aparentemente la data esta limpia, es decir, no hay errores de captura o correcciones neesarias
conteo de observaciones por tienda y por estado
dfRaw |>
count(estado, tienda)
NA
NA
Veamos como esta la informacion por estados
dfRaw |>
count(estado)
tenemos informacion de 23 estados
numero de tiendas por estados ordenados de mayor al menor
dfRaw |>
count(estado,tienda) |>
count(estado) |>
arrange(desc(n)) |>
kbl() |>
kable_classic(full_width=F, html_font = "Cambria")
| estado | n |
|---|---|
| Nuevo León | 4 |
| Puebla | 4 |
| Sinaloa | 4 |
| Veracruz de Ignacio de la Llave | 4 |
| Chiapas | 3 |
| Ciudad de México | 3 |
| Guanajuato | 3 |
| México | 3 |
| Sonora | 3 |
| Chihuahua | 2 |
| Guerrero | 2 |
| Querétaro | 2 |
| Tabasco | 2 |
| Zacatecas | 2 |
| Baja California | 1 |
| Baja California Sur | 1 |
| Coahuila de Zaragoza | 1 |
| Durango | 1 |
| Michoacán de Ocampo | 1 |
| Morelos | 1 |
| Oaxaca | 1 |
| Quintana Roo | 1 |
| Yucatán | 1 |
vemos que nuevo leon, puebla, sinaloa y veracruz son los que tienen mas tiendas, teniendo 4
ver los diferentes valores de cada variable
map(.x = dfRaw,.f = \(x) x |> as_tibble() |> distinct())
$...1
$Fecha
$Segmento
$hora_llegada
$hora_llamado
$hora_salida
$caja
$tienda
$status
$estado
NA
vemos que hay 50 tiendas, 77 cajas aunque esas cajas en realidad no son unicas 3 segmentos tenemos informacion de 396 dias
Cuales cajas hay?
dfRaw |>
count(caja)
como era de esperarse, no todas las tiendas tienen el mismo numero de cajas y en la mayoria de las tiendas se observa que hay mas cajas de retail y bancos
Conclusion
Parece que es mejor explorar la informacion dividida esta por grupos, la division natural que se me ocurre es: estado > tienda > fecha (dia)
A modo de muestra escogeremos solo el estado de puebla y lo revisaremos
dfRaw |>
filter(estado == 'Puebla') |>
group_by(tienda) |>
summarise(n())
La tienda con mayor numero de observaciones es la BA
Veamos cual de las tiendas de puebla tiene mayornumero de cajas
dfRaw |>
filter(estado == 'Puebla') |>
group_by(tienda, caja) |>
summarise(n=n()) |>
count(tienda)
`summarise()` has grouped output by 'tienda'. You can override using the `.groups` argument.
como era de esperase, la tienda mas grande ‘tienda_BA’ tiene mayor numero de cajas
Ver los tipos de cajas que existen en estas tiendas de puebla
dfRaw |>
filter(estado == 'Puebla') |>
count(caja, tienda) |>
count(caja)
NA
no todas las tiendas en puebla tienen promotoria, algunas solo cajas y ventanillas
dfRaw |>
filter(estado == 'Puebla', tienda == 'Tienda_BA') |>
count(caja)
en esta tienda solo hay 37 cajas
Veamos las observaciones de una sola tienda
dfRaw |>
filter(estado == 'Puebla', tienda == 'Tienda_BA')
vemos que hay algunas observaciones que el segmento es diferente al caja en la que fue atendido como se habia explicado en las sesiones de teams
hay unas cajas que tienen su nombre en mayusculas, causa intriga que no son tantas observaciones, veamos esas observaciones
dfRaw |>
filter(caja %in% paste0('Caja_', letters[1:5]))
dfRaw |>
filter(caja %in% paste0('Caja_', letters[1:5])) |>
count(estado)
dfRaw |>
filter(caja %in% paste0('Caja_', letters[1:5])) |>
count(tienda)
dfRaw |>
filter(estado == 'Ciudad de México') |>
distinct(caja)
dfRaw |>
filter(estado == 'Ciudad de México') |>
count(Segmento, caja) |>
View()
Conclusion:
parece que en realidad todas las observaciones pertenecen a un solo estado, ay dos observaciones que son de otros dos estados, podria deberse a un error debido a la forma en como se generaron los datos
y tambien parece que todas pertenecen a la misma tienda, supongo que fueron los datos que en la simulacion se generaron primero y luego en los siguientes se cambio la instruccion para que las cajas sean con minuscula.
parece que las observaciones las generaron por tienda por el patron que se tiene
Creo recomendable estudiar esta tienda por separado, ya que no parece seguir el mismo comportamiento de las otras tiendas Tiena_AA ciudad de mexico
Derivado de nuestro analisis exploratorio, y de las sesiones de dudas en vivo del reto coppel, concluimos despues de un cuidadoso razonamiento que el proceso fisico que se lleva a cabo en las tiendas coppel es el siguiente:
Al llegar a la tienda, dependiendo de su tramite, el cliente es mandado a una fila en donde esperara a ser atendido, esta fila puede ser alguno de los siguientes tres segmentos: banco, retail o promotoria, y aunque el cliente pueda ser dirigido a una fila, y terminar su proceso en otro, cuando un cliente es redirigido a una nueva fila, este cliente tiene que esperar a ser atendido en este nuevo segmento, (esto lo concluimos despues de analizar esas observaciones en donde el cliente inicia su proceso en un lado y termina en otro), es decir, en realidad cada uno de estos segmentos corresponde a un sistema indpendiente de lineas de espera. Y llas observaciones que inician su tramite en un lado pero terminan en otro, en realidad perteecen a este ultimo, ya que es en esa fila en donde esta esperando y en donde contribuye a que aumenten los tiempos de espera para los otros clientes.
Por tanto se procede a crear una nueva variable con el segmento real al que pertenecen, que es en base a la caja a la que realiza su tramite, caja para retail, ventanilla para banco y p para promotoria.
Asi, cada uno de estas filas son procesos independientes y pueden corresponder a un modelo de colas G/G/C
Es decir, para cada dia, en cada tienda se llevan a cabo de manera independiente tres procesos de filas de espera.
igualmete se hara una conversion a hora sexagesimal para poder entender de mejor manera los horarios, ya que aunque en el analisis no influya el como esten representadas estas horas, para el entendimiento humano si hara mas sentido si los tiempos los tenemos en formal natural.
dfRaw |>
mutate(hora_llegada_sexagesimal = as_sexagesimal(hora_llegada),
hora_llamado_sexagesimal = as_sexagesimal(hora_llamado),
hora_salida_sexagesimal = as_sexagesimal(hora_salida)) |>
ggplot(aes(x = hora_llegada_sexagesimal)) +
geom_histogram()
Parece que hay algunas observaciones con horarios de llegada fuera de un horario de atencion normal en una tienda
veamos donde estan la mayoria de las observaciones
dfRaw |>
mutate(hora_llegada_sexagesimal = as_sexagesimal(hora_llegada),
hora_llamado_sexagesimal = as_sexagesimal(hora_llamado),
hora_salida_sexagesimal = as_sexagesimal(hora_salida)) |>
mutate(hora_bin = hour(hora_llegada_sexagesimal)) |>
count(hora_bin)
dfRaw |>
slice_sample(prop = 0.01) |>
mutate(hora_llegada_sexagesimal = as_sexagesimal(hora_llegada),
hora_llamado_sexagesimal = as_sexagesimal(hora_llamado),
hora_salida_sexagesimal = as_sexagesimal(hora_salida)) |>
mutate(hora_bin = hour(hora_llegada_sexagesimal)) |>
nest(.by = c(estado,tienda)) |>
glimpse()
Rows: 50
Columns: 3
$ estado <chr> "Coahuila de Zaragoza", "Sinaloa", "Sonora", "Quint…
$ tienda <chr> "Tienda_AO", "Tienda_H", "Tienda_AZ", "Tienda_B", "…
$ data <list> [<tbl_df[4282 x 12]>], [<tbl_df[4519 x 12]>], [<tb…
Aqui estamos particionando la data de acuerdo al estado y tienda
Debido a que no cuento con una computadora muy poderosa, tomare una muestra de unicamente dos estados para hacer algunos test, pero estos analisis perfectamente se pueden hacer con todas las demas observaciones
dfSample <- dfRaw |>
filter(estado == 'Puebla' | estado == 'Sonora') |>
mutate(hora_llegada_sexagesimal = as_sexagesimal(hora_llegada),
hora_llamado_sexagesimal = as_sexagesimal(hora_llamado),
hora_salida_sexagesimal = as_sexagesimal(hora_salida)) |>
mutate(hora_llegada_bin = hour(hora_llegada_sexagesimal))
dfSample
NA
NA
nos quedaremos con las observaciones que tienen logica, escribir aqui la justificacion de porque retiramos estas IMPORTANTE
ver cuantas observaciones voy a quitar En el caso de hora llegada, la mayoria de estas observaciones hora llegada fuera de horario de tienda, tienen una hora de salida casi justo al horario de apertura de la tienda, por lo cual se puede deducir que al inicio del dia, se revisa qe no haya tickets en espera y se cancelan esos
dfSample <- dfSample |>
filter( (hora_llegada<=hora_llamado & hora_llamado<=hora_salida & hora_llegada<=hora_salida) ) |>
filter(hora_llegada_bin >= 8 & hora_llegada_bin<=22) |>
mutate(
SEGMENTO_FISICO = case_when(
str_starts(str_to_lower(caja), 'caja') ~ 'retail',
str_starts(str_to_lower(caja), 'p') ~ 'afiliacion',
str_starts(str_to_lower(caja), 'ventanilla') ~ 'banco',
.default = 'error'),
dia_semana = wday(Fecha, label = T, abbr = F),
dias_finde =
if_else(dia_semana %in% c('viernes', 'sábado', 'domingo'),
'fin de semana', 'entre semana'))|>
rename(Id=...1) |>
relocate(estado, tienda, Fecha, dia_semana, dias_finde,
hora_llegada_sexagesimal, hora_llamado_sexagesimal,
hora_salida_sexagesimal,
Segmento, SEGMENTO_FISICO, caja, status, hora_llegada,
hora_llamado, hora_salida, Id)
dfSample
Para poder analizar los datos de las llegadas de los clientes es necesario observar su comportamiento
dfSampleClean <- dfSample |>
group_by(estado, tienda, Fecha, SEGMENTO_FISICO) |>
arrange(hora_llegada, .by_group = T) |>
mutate(interarribo =
hora_llegada - lag(hora_llegada)) |>
mutate(
tiempo_espera = hora_llamado - hora_llegada,
tiempo_servicio = hora_salida - hora_llamado) |>
ungroup() |>
mutate(interarribo_sexagesimal = as_sexagesimal(interarribo),
tiempo_espera_sexagesimal = as_sexagesimal(tiempo_espera),
tiempo_servicio_sexagesimal = as_sexagesimal(tiempo_servicio)) |>
relocate(estado, tienda, Fecha, dia_semana, dias_finde, Segmento, SEGMENTO_FISICO,
hora_llegada_sexagesimal, interarribo_sexagesimal,
hora_llamado_sexagesimal, tiempo_espera_sexagesimal,
hora_salida_sexagesimal, tiempo_servicio_sexagesimal,
caja, status, hora_llegada,
hora_llamado, hora_salida, , hora_llegada_bin, Id)
glimpse(dfSampleClean)
Rows: 2,338,465
Columns: 23
$ estado <chr> "Puebla", "Puebla", "Puebla", "P…
$ tienda <chr> "Tienda_AY", "Tienda_AY", "Tiend…
$ Fecha <date> 2023-01-26, 2023-01-26, 2023-01…
$ dia_semana <ord> jueves, jueves, jueves, jueves, …
$ dias_finde <chr> "entre semana", "entre semana", …
$ Segmento <chr> "afiliacion", "afiliacion", "afi…
$ SEGMENTO_FISICO <chr> "afiliacion", "afiliacion", "afi…
$ hora_llegada_sexagesimal <time> 11:55:40, 13:38:51, 15:55:58, 1…
$ interarribo_sexagesimal <time> NA, 01:43:11, 02:17:07, 0…
$ hora_llamado_sexagesimal <time> 18:29:31, 18:19:05, 18:07:04, 1…
$ tiempo_espera_sexagesimal <time> 06:33:51, 04:40:14, 02:11:06, 0…
$ hora_salida_sexagesimal <time> 18:31:48, 18:22:31, 18:15:31, 1…
$ tiempo_servicio_sexagesimal <time> 00:02:17, 00:03:26, 00:08:27, 0…
$ caja <chr> "p_h", "p_g", "p_m", "p_m", "p_h…
$ status <chr> "Atendido", "Atendido", "Atendid…
$ hora_llegada <dbl> 0.4969954, 0.5686505, 0.6638704,…
$ hora_llamado <dbl> 0.7705023, 0.7632569, 0.7549120,…
$ hora_salida <dbl> 0.7720880, 0.7656412, 0.7607801,…
$ hora_llegada_bin <int> 11, 13, 15, 16, 18, 18, 18, 18, …
$ Id <dbl> 317281, 317285, 317260, 317292, …
$ interarribo <dbl> NA, 7.165509e-02, 9.521991e-02, …
$ tiempo_espera <dbl> 0.273506944, 0.194606481, 0.0910…
$ tiempo_servicio <dbl> 1.585648e-03, 2.384259e-03, 5.86…
dfSampleClean
dfSampleClean |>
filter(is.na(interarribo_sexagesimal))
JUSTIFICACION DE PORQUE SE DEJAN LAS OTRAS VARIABLES INTACTAS PARA LLAS SALIDAS, PUEDE QUE SEAN VERDADERAS O FALSAS, YA QUE ALGUN CLIEMTE PUDO QUEDARSE EN LA TIENDA HASTA QUE LE RESOLVIERAN SU cSITUACION Y GENUINAMENTE ESAS OBSERVACIONES REPRESENTEN A ESOS VLIENTES, SIN EMBARGO LO MAS PROBABLE ES QUE SEAN TICKETS QUE A LOS CAJEROS SE LES OLVIDO MARCAR SU HORA DE SALIDA, SIENDA HASTA HORAS DESPUES QUE SE ‘CIERRAN’ AQUELLOS TICKETS, LAMENTABLEMENTE UNICAMENTE CON LA INFORMACION DE LOS DATOS NO HAY FORMA DE DETERMINAR CUALES SI Y CUALES NO PERTENECEN A CADA UNO, SIN EMBARGO ESTAS OBSERVACIONES SI DAN INFORMACION ACERCA DE LOS TIEMPO DE LLEGADA, INTERARRIVO cDE LOS CLIETES, ASI COMO DEL TIEMPO EN ESPERA EN LA FILA, POR TANTO, ESTAS OBSERVACIONES SE DECIDE DEJARLAS TAL CUAL ESTAN
DISTRIBUCION DE LOS DATOS DE ACUERDO AL DIA DE LA SEMANA
dfSampleClean |>
ggplot(aes(x = dia_semana)) +
geom_bar()
DISTRIBUCION DE ACUERDO A LA HORA DEL DIA
JUSTIFICAR EL PORQUE DE ESTA DIVISION CON LOS ARTICULOS QUE MANDE
dfSampleClean |>
ggplot(aes(x = hora_llegada_bin)) +
geom_bar()
# veamos al final con cuanto nos quedamos por grupo
dfSampleClean |>
count(hora_llegada_bin)
dfSampleClean |>
group_by(estado, tienda, Fecha, SEGMENTO_FISICO) |>
summarise(promedio_espera = as_sexagesimal( mean(tiempo_espera)))
`summarise()` has grouped output by 'estado', 'tienda', 'Fecha'. You can override using the `.groups` argument.
Con la justificacion dada al el inicio del procesamiento de datos, se contruye la siguiente estructura de datos en la cual se anidan conjuntos de datos para cada grupo, permitiendo asi poder estudiar cada uno de estos grupos (particion) de manera simultanea
# Particion datos
dfSampleAgrupado <- dfSampleClean |>
nest(.by = c(estado, tienda, Fecha, SEGMENTO_FISICO))
dfSampleAgrupado |>
glimpse()
Rows: 7,656
Columns: 5
$ estado <chr> "Puebla", "Puebla", "Puebla", "Puebla", "P…
$ tienda <chr> "Tienda_AY", "Tienda_AY", "Tienda_AY", "Ti…
$ Fecha <date> 2023-01-26, 2023-01-26, 2023-01-26, 2023-…
$ SEGMENTO_FISICO <chr> "afiliacion", "banco", "retail", "afiliaci…
$ data <list> [<tbl_df[19 x 19]>], [<tbl_df[25 x 19]>],…
Con los datos agrupados de la forma anterior, queda comprobar si nuestra variable de interarribo y tiempo de servicio siguen una distribucion exponencial, ya que de ser asi, podemos ajustar nuestro sistema a un modelo M/M/c, donde el numero de cajas es variable para cada tienda.
ASI QUE PARA CADA GRUPO COMPROBAREMOS QUE SI EFECTIVAMENTE SE SIGUE ESA DISTRIBUCION
fun_pvalue_serv <- function(df) {
ajExp <- fitdist(as.numeric(na.omit(df$tiempo_servicio)), 'exp', method = 'mle')
test <- ks.test(as.numeric(na.omit(df$tiempo_servicio)), 'pexp', rate = ajExp$estimate)
return(test$p.value)
}
fun_rate_serv <- function(df) {
ajExp <- fitdist(as.numeric(na.omit(df$tiempo_servicio)), 'exp', method = 'mle')
return(ajExp$estimate)
}
fun_pvalue_inter <- function(df) {
ajExp <- fitdist(as.numeric(na.omit(df$interarribo)), 'exp', method = 'mle')
test <- ks.test(as.numeric(na.omit(df$interarribo)), 'pexp', rate = ajExp$estimate)
return(test$p.value)
}
fun_rate_inter <- function(df) {
ajExp <- fitdist(as.numeric(na.omit(df$interarribo)), 'exp', method = 'mle')
return(ajExp$estimate)
}
dfSampleAgrupado2 <-
dfSampleAgrupado |>
slice_sample(prop = 0.05) |>
mutate(
rate_interarribo =
map(data,
possibly(fun_rate_inter)),
rate_servicio =
map(data, possibly(fun_rate_serv)),
p_value_interarribo =
map(data,
possibly(fun_pvalue_inter)),
# p_value_servicio =
# map(data,
# possibly(fun_pvalue_serv)),
n_cajas =
map(data, possibly(\(x) {x |> distinct(caja) |> nrow() }) )
)
Aviso: There were 2625 warnings in `mutate()`.
The first warning was:
ℹ In argument: `rate_interarribo = map(data,
possibly(fun_rate_inter))`.
Caused by warning in `dexp()`:
! Se han producido NaNs
ℹ Run ]8;;ide:run:dplyr::last_dplyr_warnings()dplyr::last_dplyr_warnings()]8;; to see the 2624 remaining
warnings.
dfSampleAgrupado2 |>
glimpse()
Rows: 382
Columns: 9
$ estado <chr> "Puebla", "Sonora", "Sonora", "Puebla"…
$ tienda <chr> "Tienda_J", "Tienda_F", "Tienda_C", "T…
$ Fecha <date> 2023-12-19, 2023-07-29, 2023-06-29, 2…
$ SEGMENTO_FISICO <chr> "banco", "banco", "afiliacion", "retai…
$ data <list> [<tbl_df[616 x 19]>], [<tbl_df[348 x …
$ rate_interarribo <list> 1441.094, 844.5058, 93.93573, 1632.24…
$ rate_servicio <list> 671.1018, 502.4599, 487.7632, 433.082…
$ p_value_interarribo <list> 0.001493162, 0.3793288, 0.5187158, 3.…
$ n_cajas <list> 3, 3, 2, 6, 2, 5, 3, 8, 7, 3, 6, 1, 2…
Podemos ver que efectivamente, para cada uno de esos grupos, que corresponde a un segmento en un dia, para cada tienda, en cada estado, se sigue una distribucion exponencial de los datos de interarribo y tiempos de servicio
dfSampleAgrupado2 <- dfSampleAgrupado2 |>
unnest(n_cajas, p_value_interarribo, rate_servicio, rate_interarribo, data)
Aviso: `unnest()` has a new interface. See `?unnest` for details.
ℹ Try `df %>% unnest(c(n_cajas, p_value_interarribo, rate_servicio, rate_interarribo,
data))`, with `mutate()` if needed.
ENTONCES AJUSTAREMOS UN MODELO M/M/C
Con base en lo anterior, y en los siguietes articulos, vemos que el comportamiento de los clietes varia a lo largo del a;o, no es lo mismo los tiempos de espera en la ma;ana que en la tarde, asi como tampoco los fines de semana o entre semana.
Como el objetivo de nuestro estudio es encontrar una forma de reducir los tiempos de espera, concluimos que la mejor forma de modelar estos sistemas es analizando el comportamiento por hora del dia, y dia de la semana (mencianar articulos cientificos), asi, para cada grupo conformado por una hora del dia, un dia de la semana, para cada tienda, obtendremos su modelo
Esto lo haremos igual que usando ‘nested data frames’ y aplicando el modelo para cada grupo
df22 <- dfSampleAgrupado2 |>
group_by(estado, tienda, dia_semana, SEGMENTO_FISICO, hora_llegada_bin) |>
summarise(rate_interarribo = mean(rate_interarribo, na.rm=T, trim=0.05),
rate_servicio = mean(rate_servicio, na.rm=T, 0.05),
n_servers = round(mean(n_cajas, na.rm=T)),
median_espera = median(tiempo_espera))
`summarise()` has grouped output by 'estado', 'tienda', 'dia_semana', 'SEGMENTO_FISICO'. You can override using the `.groups` argument.
df22
Ahora ajustamos el modelo
df33 <- df22 |>
ungroup() |>
# nest(.by = c(estado, tienda, dia_semana, hora_llegada_bin)) |>
nest(rate_interarribo,rate_servicio,n_servers) |>
mutate(modelMMC= map(data, possibly(\(x) {
input <- NewInput.MMC(lambda = x$rate_interarribo, mu = x$rate_servicio, c = x$n_servers)
model <- QueueingModel(input)
return(model)
} ))) |>
mutate(Resumen_model = map(data, possibly(\(x) {
input <- NewInput.MMC(lambda = x$rate_interarribo, mu = x$rate_servicio, c = x$n_servers)
model <- QueueingModel(input)
resumen <- summary(model)
return(as_tibble(resumen$el))
} ))) |>
unnest(Resumen_model) |>
mutate(median_espera_sexagesimal = as_sexagesimal(median_espera),
Wq_sexagesimal = as_sexagesimal(Wq))
Aviso: Supplying `...` without names was deprecated in tidyr 1.0.0.
Please specify a name for each selection.
Did you want `data = c(rate_interarribo, rate_servicio, n_servers)`?
Throughput is: 181.725689688908
Utilization exceeds 100% use!!: 128.013425004069%
Throughput is: 227.646989078742
Utilization exceeds 100% use!!: 106.738154082154%
Throughput is: 299.251420100899
Utilization exceeds 100% use!!: 103.237743531182%
Throughput is: 237.651824829999
Utilization exceeds 100% use!!: 111.534174472582%
Throughput is: 181.725689688908
Utilization exceeds 100% use!!: 128.013425004069%
Throughput is: 227.646989078742
Utilization exceeds 100% use!!: 106.738154082154%
Throughput is: 299.251420100899
Utilization exceeds 100% use!!: 103.237743531182%
Throughput is: 237.651824829999
Utilization exceeds 100% use!!: 111.534174472582%
glimpse(df33)
ESTO ES LO MAS IMPORTANTE
Con esta estructura, en donde para cada hora del dia y cada dia de la semana en cada tienda tenemos el modelo que describe como se comporta el sistema de filas de espera en cada segmento, podemos entonces buscar una forma de optimizar estos tiempos,
para esto nos basaremos en la variable que creamos de tiempo de espera limite, que es el tiempo a partir del cual los cliemtes empiean a abandonar las filas
(aqui hacer una mejor explicacion)
Para lograr esto, modelaremos para cada rupo como se comportan el sistema con mas o menos servidores y observaremos su tiempo de espera en fila (Wq) (es resultado del modelo)
por ejemplo, para un caso particular que corresponde a la tienda AY del estado de puebla, en un dia domingo, a las 19:00 hrs del dia se tiene que:
i <- 10
# df33$modelMMC[[i]] |> summary()
modeltemp <- df33$modelMMC[[i]]
# modeltemp$Inputs$lambda
# modeltemp$Inputs$mu
# summary(modeltemp)
# df33 |> View()
# OPTIMIZAR UN MODELO
servers <- 1:10
ltemp <- map_vec(servers,possibly(\(x) {
inp <- NewInput.MMC(lambda = modeltemp$Inputs$lambda,
mu = modeltemp$Inputs$mu,
c = x)
model22 <- QueueingModel(inp)
return(Wq(model22))
}, otherwise = NA) )
Throughput is: 480.973034429974
Utilization exceeds 100% use!!: 314.075704596105%
Throughput is: 961.946068859949
Utilization exceeds 100% use!!: 157.037852298053%
Throughput is: 1442.91910328992
Utilization exceeds 100% use!!: 104.691901532035%
# str(ltemp)
# View(ltemp)
df55 <- tibble(
servers=servers,
Wq=ltemp
)
Tenemos aqui sus correspondientes tiempos de espera dependiendo del numero de servidores
df55 |>
kbl() |>
kable_classic(full_width = F)
| servers | Wq |
|---|---|
| 1 | NA |
| 2 | NA |
| 3 | 0.0066491 |
| 4 | 0.0006731 |
| 5 | 0.0001635 |
| 6 | 0.0000443 |
| 7 | 0.0000119 |
| 8 | 0.0000031 |
| 9 | 0.0000007 |
| 10 | 0.0000002 |
ASi podemos ahora compararlo con el tiempo de espera limite, (el cual calculamos con la mediana) y explicar el porque de la madiana.
En general no tiene que ser obligatoriamente ese criterio, puede ser cualquier otro. por eso es que se divide nuestro estudio en 5 partes
df55 |>
na.omit() |>
ggplot(aes(x = servers, y = as_sexagesimal(Wq))) +
geom_line(color='blue') +
geom_hline(yintercept = df33$median_espera_sexagesimal[i], color='red')
NA
NA
Ahora graficamos el desempe;o de estos modelos (recordar que esto es solo para uno de los grupos, osea las 7pm en domingo en la tienda AY)
Vemos que el numero de servidores actual es:
df33$c[i]
[1] 5
Sin embargo vemos que el numeor optimo de servidores es el numero maximo que este por debajo del tiempo de espera limite que en este caso es de 6 servidores.
Conclusion
Este forma de modelar los datos, nos permite ajustar para cada hora del dia el numero optimo de cajas a tener en servicio para reducir en numero de climetes que abandonan las filas.